namespace KrigingTest;

public class FastIDW
{
    public FastIDW(List<DataRecord> datas)
    {
        Datas = datas;
    }

    private double[] Lons;
    private double[] Lats;
    private double[] Values;

    public List<DataRecord> Datas;
    public int MaxEffectPoints { get; set; } = 200;

    public double[,] Interpolate2Grid(double startLon, double endLon, double startLat, double endLat,
        double maxDistance, double delta = 0.01f)
    {
        // 计算网格的行数和列数
        int cols = (int)Math.Round((endLon - startLon) / delta, 0) + 1;
        int rows = (int)Math.Round((endLat - startLat) / delta, 0) + 1;

        double[,,] wgrid = new double[rows, cols, MaxEffectPoints];
        double[,,] vgrid = new double[rows, cols, MaxEffectPoints];
        double[,] grid = new double[rows, cols];
        int[,] effectPoints = new int[rows, cols];

        // 边距设置
        double margin = 0.05; // 边距大小，可以根据需要调整

        // 创建影响权重矩阵
        double[,] effectWei = CreateEffectWeight(maxDistance, out int weihalfSize, delta);

        for (int i = 0; i < Datas.Count; i++)
        {
            // 过滤掉不在网格范围内的点
            if (Datas[i].Lon < startLon - margin || Datas[i].Lon > endLon + margin ||
                Datas[i].Lat < startLat - margin || Datas[i].Lat > endLat + margin) continue;

            // 确定站点所在的网格
            int colIndex = (int)Math.Round((Datas[i].Lon - startLon) / delta);
            int rowIndex = (int)Math.Round((Datas[i].Lat - startLat) / delta);
            if (colIndex < 0 || colIndex >= cols || rowIndex < 0 || rowIndex >= rows) continue;

            // 记录权重和站点值
            wgrid[rowIndex, colIndex, 0] = -1;
            vgrid[rowIndex, colIndex, 0] = Datas[i].WindSpeed;

            // 更新影响权重矩阵
            for (int m = -weihalfSize; m <= weihalfSize; m++)
            {
                for (int n = -weihalfSize; n <= weihalfSize; n++)
                {
                    if (m == 0 && n == 0) continue;
                    int colEff = colIndex + m;
                    int rowEff = rowIndex + n;
                    if (colEff < 0 || colEff >= cols || rowEff < 0 || rowEff >= rows) continue;

                    int level = effectPoints[rowEff, colEff];
                    if (level >= MaxEffectPoints) continue;

                    wgrid[rowEff, colEff, level] = effectWei[m + weihalfSize, n + weihalfSize];
                    vgrid[rowEff, colEff, level] = Datas[i].WindSpeed;
                    effectPoints[rowEff, colEff]++;
                }
            }
        }

        // 计算网格插值
        for (int i = 0; i < rows; i++)
        {
            for (int j = 0; j < cols; j++)
            {
                double totalW = 0;
                double wei = 0;
                bool staPoint = false;

                for (int k = 0; k < effectPoints[i, j]; k++)
                {
                    if (wgrid[i, j, k] == -1)
                    {
                        grid[i, j] = vgrid[i, j, k];
                        staPoint = true;
                        break;
                    }

                    totalW += wgrid[i, j, k] * vgrid[i, j, k];
                    wei += wgrid[i, j, k];
                }

                if (!staPoint)
                {
                    grid[i, j] = wei != 0 ? totalW / wei : 9999f;
                }
            }
        }

        return grid;
    }

    /// <summary>
    /// 创建一个站点对网格的权重矩阵
    /// </summary>
    /// <param name="distance">影响距离</param>
    ///  <param name="weiHalfSize">权重矩阵的大小</param>
    /// <param name="delta">网格间隔</param>
    private double[,] CreateEffectWeight(double distance, out int weiHalfSize, double delta = 0.01f)
    {
        weiHalfSize = (int)(distance / delta);
        int weiSize = weiHalfSize * 2 + 1;
        double[,] weight = new double[weiSize, weiSize];
        for (int i = 0; i < weiSize; i++)
        {
            for (int j = 0; j < weiSize; j++)
            {
                double dis2 = (i - weiHalfSize) * (i - weiHalfSize) + (j - weiHalfSize) * (j - weiHalfSize);

                double w = 0;
                double dis = distance / delta;
                if (dis2 > dis * dis)
                {
                    w = 0;
                    continue;
                }

                if (dis2 != 0)
                    w = 1f / Math.Sqrt(dis2);
                else
                    w = -1;
                weight[i, j] = w;
            }
        }

        return weight;
    }

    public double CalculateMaxEffectDistance(double xmin, double xmax, double ymin, double ymax, int pointCount)
    {
        if (pointCount <= 0)
        {
            throw new ArgumentException("点数必须大于0", nameof(pointCount));
        }

        double area = (xmax - xmin) * (ymax - ymin);
        double density = pointCount / area;

        // 使用密度的平方根的倒数来估算maxEffectDistance
        double estimatedDistance = 1 / Math.Sqrt(density);

        // 应用一个缩放因子来调整结果
        double scaleFactor = 1.0; // 这个值可以根据需要调整
        double result = estimatedDistance * scaleFactor;

        return Math.Max(0.01, Math.Round(result, 2));
    }
}